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1 Introduction 

An inverse problem needs to be solved in order to 
localise the brain activity that underlies the recorded 
MagnetoEncephalogram (MEG). However, this 
inverse problem is non-unique. The non-uniqueness 
can be reduced by constraining solutions to lie on 
and be orientated perpendicular to the cortical 
surface. In practice however, the accuracy to which 
the true position and orientation of the surface can be 
determined is limited by MEG-MRI co-registration 
errors. The FOCUSS [1] and constrained dipole fit 
algorithms used in this study incorporated constraints 
by means of a source space. This source space was 
constructed by placing equivalent current dipoles on 
a mesh defining a simulated cortical fold. The effects 
of the choice of the mesh-spacing and errors in the 
location and local orientation of the source space 
were examined. 

Simulated data were used so that the exact location 
of the target sources was known and the quality of 
the source reconstructions could therefore be 
objectively evaluated. 

2 Methods 

2.1 Source space 

The 3D source space was constructed by 
concatenating a set of parallel contours. The source 
space consisted of two parallel sulci, 4cm deep, 
joined by a gyrus, approximating a true cortical fold 
(Fig. 1). A source space very similar in overall shape 
was used in a simulation study by Genger and 
Williamson [2], which was an extension of the 
model of a single sulcus used by Wang [3]. 
Mesh-spacing was varied by downsampling an 
original mesh (mesh-spacing of 0.5mm). The source 
space was formed by placing an equivalent current 
dipole at each node of the mesh. The orientation of 
the dipoles was constrained by the orientation of the 
surface. 

The source space was positioned in a single sphere 
and the configuration of a 151 channel CTF Omega 
system (CTF Systems, Inc.) was used as a simulated 
sensor array (Fig. 1). 


2.2 Target source 

A target source with a Gaussian amplitude 
distribution (g= 10mm) was centered on the gyrus in 
the original mesh (mesh-spacing 0.5mm). Modelling 
such an extended source as a single dipole should not 
introduce large errors in the source reconstructions 
[4]. The source amplitude was sinusoidally 
modulated in time so that one cycle fitted in the 
selected time-window (50 samples), 
a) b) 



Figure 1: a) Source space modeled as a cortical fold 
viewed from the top (a) and right-back (b). The 
channel array is depicted by the green filled circles 
and the center of the target source on the gyrus is 
indicated by an arrow. A= Anterior; P=Posterior 

2.3 Noise 

Zero mean Gaussian noise (g= 20fT) was added to 
the simulated data in order to probe the limits of the 
tested algorithms. 

The signal-to-noise ratio (SNR) was varied by 
varying the amplitude of the target source, and hence 
the amplitude of the simulated MEG data. 

Monte Carlo simulations were performed in order to 
assess the stability of the source reconstructions. It 
was found in initial simulations that the mean 
distance from the solutions to the target source and 
the standard deviation of this distance stabilised 
when the number of realisations was in the order of 
30. 

Each dataset, containing a different noise realisation, 
was used for the reconstruction of the underlying 
sources with three distinct algorithms. The 




implementation of these algorithms is described 
below. 

2.4 Implementation of source reconstruction 
algorithms 

Unconstrained dipole fitting 

A fixed dipole model was used to explain the 
simulated data, i.e. the same dipole location and 
orientation was selected for all latencies. The 
timecourse of activation was obtained from the data 
by linear estimation once the location and orientation 
of the dipole was known. 

Initial Guess A single dipole was used in the 
dipole fit algorithm. The location of the initial guess 
was randomly chosen within a sphere with a radius 
of 2cm around the target source. The initial guess for 
the orientation of the source was randomly selected. 

Cost-function The reduced chi-square measure was 
used as the cost-function. Its value was estimated for 
each latency separately and the mean over all 
latencies was used as the value of the cost-function. 
Dipoles with a magnitude larger than lOOnAm or 
that moved outside the sphere were penalised by 
squaring the value of the cost-function. 

Optimisation algorithm A multi-start optimisation 
was used in order to avoid local minima [5]. With 
each dataset, the Simplex algorithm [6] was started 
with 10 initial guesses [7]. The solution that resulted 
in the lowest value of the cost-function was chosen 
as the solution for the dataset, referred to as the 
unconstrained dipole solution (UCDIP). 

Constrained dipole fitting 

Again a fixed dipole model was used to explain the 
data. However, this time the location of the dipole 
was restricted to lie within a distance equal to the 
mesh-spacing from the nodes in the source space. 
Sources that moved further away from the surface 
were penalised by squaring the value of the cost- 
function. Additionally, the orientation of a dipole 
was constrained to be the same as the orientation of 
the normal of the nearest node in the source space. 

Apart from these constraints, the same 
implementation as for the unconstrained dipole 
fitting algorithm was used. For each dataset this 
resulted again in 10 solutions due to the multi-start 
approach, from which the one with the lowest value 
for the cost-function was chosen, referred to as 
CDIP. Also, the nearest node in the mesh from CDIP 
was stored as CDIP Node. The latter solution was 
intended to be a fair comparison for solutions 


obtained with FOCUS S, since FOCUS S also 
incorporated the spatial and orientation constraints. 

Multi-latency FOCUSS algorithm 

The FOCUSS algorithm [1] was modified so that 
temporal information could be incorporated [8]. The 
FOCUSS algorithm progresses by using the solution 
from previous iterations as the weights in the current 
iteration, until the algorithm converges to a final 
solution. The algorithm can easily be adjusted for 
multi-latency data [8]; the value for each element on 
the diagonal in the weighting matrix is formed from 
the weighted contribution of the activity at all 
latencies for that element in the previous iteration. 
The weightings can be adjusted to alter the 
contribution from different latencies, based for 
example on the SNR at each latency. In this study all 
latencies were given equal weighting with the 
compound version of FOCUSS. 

Bias matrix The initial weighting matrix was 
focussed on a target area by selecting the values of 
the weights on the basis of the distance from the 
element in the source space associated with the 
element in the weighting matrix to the target source. 
The value of the weights fell off with the distance 
from the target source, following a 3D-Gaussian 
distribution (g= 10mm). This initial weighting 
matrix was then multiplied by the squared average 
distance from the elements to the sensors and by an 
additional factor to counteract the bias towards the 
sensors. This additional factor was equal to the 
inverse of the difference between the distance to the 
farthest sensor and the nearest sensor. 

Noise regularisation Noise regularisation was 
performed by means of Truncated Singular Value 
Decomposition. For each dataset, numerous 
truncation levels were explored. With each selected 
truncation level a source reconstruction was 
obtained. The final solution for each dataset was 
chosen on the basis of the lowest mean distance from 
the reconstructed source(s) to the target source. Note 
that this would not be possible with recorded data. 

2.5 Simulations 

Influence of mesh-spacing 

The spacing of the meshes was 4, 6, 8 or 10mm. 
These meshes were used with the constrained dipole 
fit algorithms for the reconstruction of the sources 
underlying the simulated data. Additionally, a mesh 
with a spacing of 2mm was used for the constrained 
dipole fits. Datasets with mean SNRs of 1.2, 1.5 and 
3 were used in the Monte Carlo simulations. 



Influence of displacement errors 

The mesh and sphere origin were shifted in a random 
direction by a distance taken from a Gaussian 
distribution (g= 2, 4, 6, 8 and 10mm) with zero 
mean. These shifts in the mesh location and the 
origin represented errors in the location of the mesh 
due to errors in the co-registration of the MEG and 
MRI co-ordinate systems. The shifted mesh was 
subsequently used by the constrained source 
reconstruction algorithms to find the solutions for 
simulated data with a mean SNR of 1.5. The spacing 
of the mesh was 6mm. 

Influence of orientation errors 

The normal at each node was perturbed by a random 
angle taken from a Gaussian distribution (g=10, 20, 
30 and 40°) with zero mean. Sources were 
reconstructed from simulated data with a mean SNR 
of 1.5, using a mesh with a spacing of 6mm. 

3 Results 

3.1 Influence of mesh-spacing 

The influence of the mesh-spacing on the solutions 
obtained with FOCUSS is illustrated in figure 2 and 
for the constrained dipole fit algorithms in figure 3. 
These figures show that: 

• For FOCUS S the mean distance error and the 
standard error increased with an increase in 
mesh-spacing. 

• The constrained dipole solutions did not depend 
on the mesh-spacing and were comparable to the 
unconstrained solutions. 

• The constrained dipole algorithm, with the final 
solution constrained to a node in the mesh, gave 
more accurate and stable solutions than the 
unconstrained dipole fit algorithm. However, the 
advantage of using constraints disappeared for 
large mesh-spacings. 

3.2 Influence of displacement errors 

Figure 4 shows that the error in the source 
reconstructions increased with increasing 
displacement errors. However, the gradient of the 
curves show that the increase in localisation error 
was only about half the increase in the displacement 
error (a 4mm displacement error resulted in an 
increase in localisation error of about 2mm). 

3.3 Influence of orientation errors 

Figure 5 demonstrates that: 

• For FOCUSS the influence of orientation errors 
on the mean localisation error was small, 


although the stability of the solution decreased 
with increasing orientation errors. 

• The accuracy of the solutions obtained with the 
constrained dipole fit algorithms decreased with 
an increase in the orientation errors. The 
decrease in accuracy was however smaller then 
3mm for orientation errors as large as 40°. Also, 
the constrained dipole fit algorithms were robust 
for orientation errors up to 10°. 



Figure 2: Influence of mesh-spacing on the mean 
localisation error and standard error for solutions 
obtained with FOCUSS, using data with different 
mean SNRs. The results for the unconstrained dipole 
fits (UCDIP) are shown for comparison. 



Figure 3: Influence of mesh-spacing on the mean 
localisation error and standard error for solutions 
obtained with constrained dipole fitting (CDIP), and 
the constrained dipole algorithm, with the final 
solution constrained to a node in the mesh 
(CDIP Node). The results for the unconstrained 
dipole fits (UCDIP) are shown for comparison. 
Data with several different mean SNRs were used. 
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Discussion 



Figure 4: Influence of displacement errors on the 
mean localisation error and standard error for 
solutions obtained with FOCUSS and the 
constrained dipole fitting algorithm (CDIP). 
CDIPNode are the solutions for the constrained 
dipole fitting algorithm, with the final solution 
constrained to a node in the mesh. The results for 
the unconstrained dipole fits (UCDIP) are shown for 
comparison. 



Figure 5: Influence of orientation errors on the 
mean localisation error and standard error for 
solutions obtained with FOCUSS and the 
constrained dipole fitting algorithm (CDIP). 
CDIP Node are the solutions for the constrained 
dipole fitting algorithm, with the final solution 
constrained to a node in the mesh. The results for 
the unconstrained dipole fits (UCDIP) are shown for 
comparison. 


• The use of anatomical constraints with a 
constrained dipole fit algorithm is only 
beneficial when the mesh-spacing is small (in 
our simulations < 4mm). 

• The increase in localisation error was only about 
half the increase in the displacement error. 

• The influence of orientation errors on the 

performance of the constrained source 

localisation algorithms was relatively small. 

Consequently, MEG-MRI co-registration errors and 
errors in the extraction of the cortical surface from a 
MRI do not have disastrous effects on the 
performance of constrained source localisation 
algorithms. 
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